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Abstract 

We calculate the subjet rates for jets produced in hadron collisions. The /cj_ algorithm 
is used to define the jets and allows the theoretical calculation to sum both the leading 
and next-to- leading logarithms in the resolution variable, y C ut- We also ensure that our 
calculation matches exactly the leading order in a s result and has sensible behaviour 
near thresholds. 
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1 Introduction 



Hadron collisions are a prolific source of hadronic jets. When these are hard and well-separated 
they can be used for a variety of high-precision tests of perturbative QCD [p]]. However, 
much of the theoretical interest comes from the soft region in which the confinement of the 
produced quarks and gluons into the observed hadrons plays a crucial role. Since this is 
outside perturbative control a number of models have been proposed to describe it. Although 
these all give reasonable descriptions of hadron-level data, to gain a deeper understanding 
one would like to bridge the gulf between the individual hadrons and the hard well-separated 
jets, in such a way that one can choose the dominant domain, smoothly moving between the 
two extremes. A powerful way of doing this is to define subjets within the jets [0. 

Subjets are defined by taking all the particles that ended up in a given jet and running a jet 
algorithm on those particles with a dimensionless resolution scale, y CU f For large y cu t, every 
jet consists of only a single subjet. As y cnt is decreased the fraction of jets that are resolved 
into two or more subjets increases, until for small y cut every jet consists of many subjets. 
Eventually for small enough y cnt , every final-state particle is resolved from every other, and 
subjet distributions become identical to hadron distributions. The most interesting region is 
the intermediate one in which the resolution scale is large enough for perturbation theory to 
be valid, but small enough for the typical multiplicity of subjets to be large. Careful study of 
this region already shows many of the features of the hadronic final state, despite being fully 
calculable in perturbation theory. 

Subjet studies are already common in e + e~ annihilation^] ]2j, |3|, f|]. However in hadron 
collisions there have been relatively few studies, mainly because most experiments use cone- 
based jet algorithms, which are not amenable to the all-orders calculations that are needed 
to make subjets interesting. However, with the advent of the k± algorithm for hadron col- 
lisions subjet studies became possible 0. k± algorithms are of the clustering type 
and are defined by a closeness measure, which decides whether two particles are clustered 
together, and a recombination scheme, which specifies how they are clustered together. Once 
these two details have been defined, the algorithm can be applied iteratively to final states 
containing any number or type (partons/hadrons/calorimeter cells) of particle. For collisions 
involving incoming hadron beams, one also has to define a closeness to the beam direction 
to ensure that the resulting cross sections obey the factorization theorem ||. The variant we 
use is the 'inclusive k± algorithm' with the 'p r recombination scheme'. 

To be precise, we define jets in the following way, in terms of a dimensionless parameter 
R, which we will see plays a radius-like role in defining the angular extent of the jets, and 

^Note that the nomenclature we use is that of hadron collisions in which one imagines a two-step procedure: 
first defining jets, then studying subjets within them. In e + e~ annihilation, this nomenclature is sometimes 
used, for example in studying the subjet structure of three-jet events, but since e + e~ annihilation is a point- 
like source of two-jet events, we can consider every jet study there to be a subjet study, since the first step is 
not needed. 
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which we assume is of order 1. We start with a list of 'particle' momenta that could consist of 
partons, hadrons or calorimeter cells, and an empty list of jet momenta. In the p t scheme, all 
particles are treated as if they were massless, so momenta are defined by only three variables, 
transverse momentum^], pt, pseudorapidity, 77, and azimuth, 0. Since these are invariant under 
Lorentz boosts along the beam direction, the whole algorithm is longitudinally boost invariant. 

1. For every pair of particles, calculate their closeness 

dij = min^py AR%, (AR% = ( Vi - Vj ) 2 + (0, - 0,) 2 ) . (1) 

Note that in the small angle limit, this reduces to the transverse momentum of the softer 
particle relative to the axis defined by the harder: 

dij « k±l = min(£ 2 , Ej) sin 2 A%, AR% < 1. (2) 

2. For every particle, i, calculate its closeness to the beam directions 

d l= pl t R 2 . (3) 

3. Find the smallest closeness <i min = min {d^, di}. 

(a) If d min = d^, particles % and j are merged, by replacing their entries in the particle 
list by a new entry, with momentum 

Pt,ij = Pt,i + Pt,j, (4) 
Vij = (pt,iVi+Pt,jVj)/Pt,ij, (5) 
<Pij = (pt,i<f>i+Pt,j(f)j)/Pt,ij- (6) 

(b) If d m [ n = di, i is said to be a completed jet: its momentum is moved from the list 
of particles to the list of jets. 

4. Continue from step 1 until the list of particles is empty. 

Note that there is an unambiguous assignment of every final-state particle to exactly one 
completed jet. Note also that although the particles are merged in order of increasing relative 
transverse momentum, it is actually their angular separation that controls whether or not 
they are merged: every merged pair within a jet has ARij < R and every jet is separated 
from every other by ARij > R. 

Having selected a particular jet k for further study, we define subjets within it, in terms 
of a dimensionless parameter y cut , by rerunning the above algorithm, but only putting those 

■'■Note that the transverse momentum, p t , is often replaced by the transverse energy, Et- Since the pt- 
scheme uses massless kinematics, the two are identical. 
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particles that were assigned to our jet in the initial list, and stopping as soon assail closenesses 
satisfy 

dij > Vcut vlk- (7) 

All particles in the list are then called subjets. 

Experimental implementations of this algorithm (see for example ||) have tended to also 
define a preclustering stage in which the calorimeter cells are clustered together with a fixed 
angular cut. We do not consider the effect of such a procedure here, so cannot directly 
compare with their data at present. In ||, and di are also defined slightly differently: 
both are divided by a factor of R 2 relative to ours. y cut is still defined through ([?]) though, 
so that it is different from ours by a factor of R 2 . This does not affect any physical results, 
except that their y cnt is therefore identical to our Y, defined in (|10|). 

In the mean number of subjets in a hadron collider jet was calculated, to leading order 
in a s and summing leading and next-to-leading logarithms of y cn t to all orders in a s . In this 
paper we go one step further by calculating the n-subjet rates, P n (y C ut), i-e. the fraction of 
jets that consist of n subjets at a given resolution scale, to the same accuracy. 

In general, one finds that P n ~ a™ -1 and that at small y cut the leading term at each order 
af, m > n— 1, is ln 2m y cut . In the region we are interested in a s In 2 y cut > 1 and order-by-order 
perturbation theory does not converge, so these (leading-logarithmic) terms must be summed 
to all orders in a s . Since each higher order in a s contributes two additional logarithms, 
we must also sum the next-to-leading logarithms, a™ ln 2m_1 y cut , before we obtain a result 
that is of leading-order accuracy (i.e. in which neglected terms are down by a power of a s , 
without any logarithmic enhancement, relative to included terms). The modified leading 
logarithmic approximation (MLLA) provides a framework to perform this next-to-leading 
logarithmic resummation (see for example p0| , [TT|j). It turns out that although the leading 
logarithms are identical in hadron collisions to the well-studied case of e + e~ annihilation, the 
next-to-leading logarithms contain essential contributions from soft gluons that are radiated 
off the incoming partons. The probability of this depends on the jet kinematics, the parton 
distribution functions, the collision type and energy and so forth, so cannot easily be treated 
analytically. However, as pointed out in 0, it is possible to manipulate all these effects into 
a single factor, which can be calculated numerically, multiplying a tower of logarithms that 
can be summed analytically. 

Finally, in order to describe the region of large y cnt , it is necessary to include in P n the 
exact tree-level contributions from the n + 1-parton final state. Although these are available 

§Note that it is possible to reach a situation in which, although (0) is satisfied, further merging would 
result in a configuration in which it was no longer satisfied, and only even later in the merging history would 
it again be satisfied. Such 'non-monoticity' is a property of the recombination scheme, and is shared by the 
only other scheme in common use, the .E-scheme. It is possible to define more complicated schemes that are 
guaranteed to be monotonic, like the -scheme, and the results in jfj] show that, as one might expect, they 
suffer smaller hadronization corrections. 
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for n up to 6 
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, we have only included them up to n = 2. Formally, our accuracy is therefore 



(1-Pi) = (l-P l )(l + 0(a s )) Wy cnt , 

P 2 = P 2 {l + 0{a s )) Vy cut , (8) 

Pn>2 = Pn>2 (1 + 0{a 8 )) y cut < 1. 

Since the higher subjet rates are small for large y C ut, this is not a severe limitation. 

Within the MLLA the flavour of a jet is well-defined, because all emission is either of a soft 
gluon, or a collinear parton pair, both of which preserve the jet's net flavour. However, after 
including the exact 0(a s ) contribution this is no longer the case as, for example, configurations 
arise in which the two subjets of a jet are both quarks. As shown in [|5J this is a tiny effect 
and for all practical purposes quark and gluon jets can be separated. Thus for all our results 
we show three distributions: for the sum of all jet types (which include the unclassifiable 
jets), and for quark and gluon jets separately. In fact one of our main results is the fact that 
the properties of quark and gluon jets of a given transverse momentum are almost entirely 
insensitive to how they were created: their position in the detector, the beam particles and 
energies, etcQ. 



2 Summing the leading and next-to-leading logarithms 

As mentioned in the introduction, the leading logarithms are identical in hadron collisions to 
in e + e~ annihilation, but the next-to-leading logarithms receive additional contributions from 



initial-state radiation. We therefore first consider the final-state contributions in Section ETT 



before adding in the initial-state contributions in Section |2.3|. This in turn leads us to a 



natural way of matching with the exact 0(a s ) matrix-element result, which we discuss in 



Section |2.4j . All of these results can be improved by an improved treatment of the threshold 



region in the resummed calculation, which we discuss in Section |2[ 



2.1 Final-state logarithms 

Within the MLLA, the hardness of a jet is characterised by the product of its energy and 
the maximum allowed angle of emission from it, in our case Q = ptR- Its subjet structure 
is resolved with respect to a cutoff on transverse momentum Q = p t \/y m t- We seek to sum 
the leading and next-to- leading logarithms of Q/Qo to all orders. To this accuracy, we are 

^This is no longer true, however, if one breaks the inclusive nature of the jet definition, by imposing cuts on 
the other jets in the event. Thus although we would expect our results to describe inclusive photoproduction 
data such as Jl3[ reasonably well, we would expect them to break down once an x-y cut is imposed, as it often 
is in experimental analyses. We consider this, and other contributions unique to photoproduction, in a future 
publication. 
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insensitive to the precise values of Q and Q provided Q/Qo is left unchanged, i.e. we are free 
to multiply both by any number of order unity. It will turn out to be useful to have Q = p t , 
and hence we define 

Q = Pt, (9) 
Qo = Pt VVcut/R = Pt VY- (10) 

Our procedure correctly sums logarithms of Y to all orders in a s , but does not keep track of 
the additional logarithms that would arise if R was much less than 1. 

Considering only final-state logarithms, the probability of finding n subjets within a jet 
of flavour a, P%(Q, Y), does not depend on how the jet was produced. In order to compute it 
we shall use a generating function, $ («, Q), defined by 



n=l 

from which it follows that 



<5> a (u,Q) = Y,u n P:(Q,Y), (11) 

n=l 

1 d n $ a (u,Q) 



p:(q,y) 



(12) 

u=0 



n\ du n 

For brevity, we suppress the dependence of the generating function on Y . The use of generat- 
ing functions in describing quark and gluon jets in e + e~ collisions is now textbook material 
see for example [JT0|, |Tl|. They obey the coupled evolution equations [|, [10], |TT 



^^hql Q' 2 Jqo/q' 

a s (z(^-z)Q') p ^ Uz) ^ fc(U) zQl) $c(U; (1 _ z)Ql) _ $a(U) Q ,fj ^ (13) 

where P a ^t> c (z) is the leading order DGLAP splitting kernel for the branching a —>■ b + c, with 
b taking an energy fraction z. In the final bracket of (|13D , the first term comes from the real 
splitting a — > b + c, while the second term comes from the corresponding virtual terms, such 
that unitarity is preserved in the sum of the two. 

These equations can be rewritten in a form in which they are easier to solve, by introducing 
the Sudakov form factors, 

a f 1 f Q " dQ' 2 fi-Qo/Q' , a B (z(l - z)Q') n , A 

A a (Q) = eX p -- £ / 2 7^ / /rv dz sV V ^> P a ^ hc {z) , (14) 
{ 2 ^ JaqI Q' a Jqo/q' 2tt J 

which sum the virtual terms to all orders, to give 

' 1 rQ 2 dQ' 2 fi-Qo/Q' 



* tt (u,Q) = M A a (Q)exp -W a %- I dz 

{ 2 Y J *Q 2 o Q 2 JQo/Q' 



g s (z(l - z)Q') „ ^ b (u,zQ')<S> c (u,(l-z)Q') \ n 

Pa^bc(z) ~ . . >. 15 

2tt $ a {u,Q') J 
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Note that the integration has soft singularities at z = and z — 1. It is easier to solve (|15|) 
if we rewrite it in such a way that all the singularities appear at z — 0. We multiply the 
integrand by 1 = z + (1 — z) and use the fact that P a ^ bc (z) = P a ^ cb (l — z) to obtain 



<& a {u,Q) = MA a (Q)expiV/ 

I 6 u 

a s {z{l - z)Q') 
2tt 



Q 2 dQ' 2 i-i-Qo/Q' 



(/2 



1 - z)P a ^ bc (z) 



$ h {u,zQ') <5> c (u,(l-z)Q') 

®a(u,Q') 



(16) 



This equation only has a soft singularity at z = 0. We are therefore free to replace (1— z) by 
1 in all smoothly- varying functions, and to rewrite the upper limit of the z integration as 1. 
The form of the soft singularity is universal so, introducing C a=q = Cf and C a=g = Ca, we 
can extract it: 



<$> a {u,Q) = MA a (Q)expj^ 
(2C a 



Q 2 dQ' 2 /-i , a s {zQ') 
dz 



I Q' 2 JQo/Q 



2tt 



%(u,zQ')+J2 



[1 - z)P a ^ bc (z) - 5 t 



2C a 



bg- 



$ b (u,zQ')$ c (u,Q') 



■ (17) 



Since we have removed the soft singularity, we can replace z by 1 in the final term. Rewriting 
somewhat, we therefore obtain: 



®a(u,Q) = uA a (Q)exp<J 



dQ' 2 a 8 (Q') 



Q 2 
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$a(u,Q') 



[ dz ({l - z)P a ^ bc (z) - dhg^^ 
Jqo/Q' V z 



Finally, the lower limit of the final z integration can be replaced by 0. Writing these 
equations explicitly, we therefore have 



$q(u,Q) = uA q (Q)exp 
® 9 {u,Q) = uA g {Q)exp 



f Q o dQ 12 r q (Q,Q')%(u,Q') 



and 



A g (Q) = exp 
A q (Q) = exp 



dQ' 2 {T 9 (Q, Q')%(u, Q>) + T f (Q') QO } 



[ Q dQ' 2 (Tg(Q,Q') + T f (Q')) 



(19) 



(20) 



Ql 
Q 2 

Ql 



dQ' 2 T q (Q,Q') 



(21) 
(22) 
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where we have defined 



17(G) 

r g (Q,Q') 
r g (g,g') 



N f a s {Q 2 ) 
6tt Q 2 
C A a s (Q' 2 ) 
2tt Q' 2 

2tt Q 



rl 



Hi'. 
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IT 

I) 



(23) 
(24) 
(25) 



Throughout this paper, we use the one-loop renormalization group equation for the running 
of a s within a given jet, but with the starting value fixed by the input parton distribution 
function set (two- loop) at scale Q — p t : 



a s (g 2 



1 + 



47T 



A) In 



Q 2 ' 
Pt 



(26) 



where & = (HCU - 2JV»/3. 

The most convenient way to solve (|19],^0|) is by power series expansion in u. Writing 



we find 



*.(«,q) = 5> b *£ b) (G). 



n=l 



*S 1} (g) 
^ 1} (g) 



A ? (g), 
A g (g), 



(27) 



(28) 
(29) 



^ 2) (g) 
*f (g) 



A 9 (g) r rfg /2 ryg,g')A 9 (g') 



A g (Q) / n rfg' 2 



r g (Q,Q')A g (Q') + r f (Q 



„ A q (Q>) 

A s (go J 



(30) 
(31) 



$f(g) 
<^ 3) (g) 



^ 2) (g) 2 

2A g (g) 



+ A 9 (g) /fdg' 2 r 3 (g,g')$f (g') 

J On 



2A g (Q) 



+ r/(g 



+ A 9 (g) | q2 dg' 2 |r 9 (g,g')$( 2) (g') 

2$? ) (g / )-§j|H$f(g / 



^ A g (gp 
A 9 (go 



(32) 



(33) 



^ - -IS + a " (q) £ dQr2 r » w - Q ' w °' m + • <34) 



& 2) (Q) 3 rQ 2 



3A g (Q) 2 ' J ^'J Q i 




+ am ■ (35) 

These results sum the large leading (a™ ln 2m Y) and next-to-leading (a™ ln 2m_1 Y) logarithms 
to all orders in a s . However they are not guaranteed to be reliable, or even physically-behaved, 
in the threshold region Y ~ 1 where logarithms of Y do not dominate. We discuss an improved 
treatment of this region in the next section. 

Having these generating functions at our disposal it is straightforward to calculate the 
n-subjet rates in hadron collisions in the 'final-state approximation' (we discuss corrections 
arising from initial state radiation in Section |2.3| ). The jet generating function $j e t is formed 
from the following admixture of the two parton types: 

$ M (U,Q) =J2 F a®a(u,Q), (36) 
a 

where the factor F a labels the fraction of events that lead to the production of a leading 
parton of type a (a = q or g), i.e. 



SXL c drfj i f 1 {x 1 )f 2 (x 2 )\M a 2 \ 2 



r °- ~ T»7m ax J„/ 1 £ /„ \ M M A A 19" \°' ) 



Here fi(xi) denotes the parton distribution function of hadron i (the factorization scale depen- 
dence and flavour dependence are both suppressed), rj is the rapidity of the recoiling parton; 
it is integrated over all phase space, i.e. r] max = ha.(y/s/p t ). The Mandelstam s is that of the 
hadron-hadron system and s is that of the parton-initiated subprocess, i.e. s = X\X 2 s. The 
2^2 matrix elements are computed to lowest order in a s . 



2.2 Threshold region 

The solutions to (|l^,^), i.e. equations (|28H35|), obey the physical boundary condition that at 
Y = 1 every jet consists of 1 subjet. However, reducing Y, one finds that $^ becomes greater 
than one, while $£ 2 - ) becomes negative, which is clearly unphysical. Only at considerably 
smaller Y do they start to become physically behaved. Formally this is not a problem, as 
our results are only strictly valid for very small Y. To get a good description of the large Y 
region one needs to match with the exact order-by-order perturbative results. 
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Phenomenologically however, it is highly desirable for the resummed results to be physi- 
cally behaved for all Y and to have thresholds at the correct point in Y. That way we expect 
the matching with fixed-order perturbation theory to be smoother, and to obtain reliable 
results at a lower perturbative order. 

There are two issues to address: the position of the threshold and the behaviour close to 
it. We tackle them in reverse order. 

On critical examination of the sequence of approximations we made to solve (|T5|), we 
find that all obey positive-definiteness of the probability distributions except the very last: 
replacing the lower limit of the z integration of the non-soft splitting functions by zero. Up 
to and including fll~8l) , the soft and non-soft parts of the splitting functions are treated on an 
equal footing, but in getting from (|l"8l) to (|l9ip0|) , we apply looser phase space limits to the 
non-soft terms than to the soft terms. Since the non-soft splitting functions are not positive 
definite, it is not surprising that extending the range over which they are integrated, without 
simultaneously increasing the range over which the soft ones are integrated, leads to negative 
distributions. 

A straightforward solution to this unphysical behaviour was proposed^ in Q : one solves 
(|j~8D without approximating the z limit. That is, one simply replaces the T functions by 



r/(Q;Q ) 
r g (Q,Q';Q ) 
T q (Q,Q';Q ) 



N f a s {Q 
2tt Q2 
C A a s {Q 12 
2tt Q' 2 
C F a s (Q' 2 
2tt Q 




T- 4 (») + 3 (t) 2 -t(f)S(f)' 




and all other results, (plp2[ , p8| - |35|) , remain unchanged. 

We now have n-subjet rates that go smoothly to zero as Y — > 1 for all n > 1 and smoothly 
to one for n — 1. However, the physical threshold for n massless partons to all be resolved 
from each other does not lie at Y — 1 but rather at Y = Y n < 1. We have not found a general 
form for Y n , but we explicitly find, for n < 4, Y n = 1/n 2 , and Y 5 = 4/(5 + 2 v /5)/5 2 « 0.422/5 2 . 
It is worth noting that even exact numerical solution of the original equation flT3| ) would not 
get the threshold position right in general. Fortuitously it does in our case for n = 2, but it 
does not, for example, in e + e~ annihilation, or in either case for larger n. 

One, completely arbitrary, way to ensure that the n-subjet rate goes smoothly to zero 
as Y — > Y n is to note that with logarithmic accuracy, Y can be replaced by any arbitrary 

fc (z) they have 



jnon- soft ( 
*bc 



II Our solution is slightly different from theirs: where we have Jg ( ^,dzP^ 

/q~/2Q' 2Q ' dz P^ b - soit (z). Although their form also gives a physically-behaved threshold at Y = I, it still 
does not treat the soft and non-soft terms equivalently. As a result it is not as good an approximation of the 
full result. 
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rescaling Y(l + 0(Y))WM. For example, if we replace Y by Y(l + 1 Y % n Y) then the n-subjet 
rates will go smoothly to zero at Y — Y n . In particular, we replace 

Pl (Y) —> P x (Y(l + 12Y)) , (41) 
P n {Y) — > P rt (Y(l + n V ~ , 2<n<4. (42) 



Canonical Jets 




ycut 



Canonical Jets 




y^cut 



Figure 1: Quark subjet rates in the NLLA 
(dotted), with the inclusion of Qo-suppressed 
terms (dashed) and with shifted thresholds 
(solid). Also shown is the sum of all subjet 
rates (n < 4). Final state radiation only. 



Figure 2: Gluon subjet rates in the NLLA 
(dotted), with the inclusion of Qo-suppressed 
terms (dashed) and with shifted thresholds 
(solid). Also shown is the sum of all subjet 
rates (n < 4). Final state radiation only. 



We show results in Figures [l] and |2|. They demonstrate that the threshold modifications 
make a considerable difference to the n-subjet rates over quite a wide range in y cut . 

Throughout this paper our canonical jets are defined to be produced in pp collisions 
with p t = 60 GeV, r] = 0, R = 1 and = 1800 GeV, and we use the CTEQ4M parton 
distribution functions [15| as implemented in PDFLIB (NSET = 34) [l6|]. In fact in the final- 
state approximation we are discussing here, the properties of quark and gluon jets depend 
only on their p t and the other parameters are relevant only for determining the quark-to-gluon 
mix. 

We will see later that, after matching with the leading order in a s results, the modifications 
make less difference, although they are certainly still not negligible. The suppression of the 
curves corresponding to the sum of the first four subjet rates at large y cut after including the 
threshold rescaling is indicative of the potentially large corrections in the threshold region. 
We shall later see that this suppression disappears after matching with the leading order in 
a s results. 
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2.3 Initial-state logarithms 



An additional complication occurs in hadron collisions, because at next-to-leading logarithmic 
accuracy gluons emitted by the incoming partons can contribute to the number of subjets 
in a jet. Since such emission must be non-collinear one logarithm is lost, so to our required 
accuracy it is sufficient to consider one initial-state gluon emission, followed by the leading- 
logarithmic evolution of that gluon, in which all subsequent emission is both soft and collinear 
to it. 

We assume that the leading-order prediction for soft initial-state gluon emission into a jet 
of flavour a is known and given by 

dP a = A a ^a s . (43) 

K ± 

The constant A a depends on the jet kinematics, the hadron collision energy and type, the 
parton distribution functions, etc. We return to how it is calculated in the next section. The 
jet's generating function, fl3"6j), then becomes 



F a $ a (w, Q) + / -±a 8 {k 2 ± )A a $ 8 (u, Q)($ g (u, k ± ) - 1) 

J Of, K i 



(44) 



The first of the additional terms (the <& 9 {u, k±) term) accounts for the evolution of the soft 
initial-state gluon from the real emission. The second (—1) term accounts for virtual correc- 
tions to the 2 — ^ 2 matrix elements and is needed in order to conserve probability. 

Since the non-collinear emission loses us one logarithm, we can work in the LLA for the 
generating functions which appear in the second and third terms on the right-hand-side of 
(f44|). The evolution equation for the generating function of a parton can then simply be 
obtained by truncating ( pf) at LLA: 

d$ a (u,Q) [Q 2 dk 2 ± a s (k 



d\nY 

Hence we can rewrite (Ell): 



rQ 2 dk 2 rv (k 2 \ 

_ C J ^±^±l^ Q){%{ ^ k±) _ iy (45) 



In writing fl4~6] ) we have avoided the need to perform the integral over the transverse mo- 
mentum of the soft gluon. We just need to compute A q and A g . Since doing this involves 
using the exact three-parton matrix-element calculation, it gives us a very convenient way of 
matching with the exact 0(a s ) result. 
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2.4 Matching with exact 0(ct s ) result 

We begin by writing down the exact matrix-element result, before going on to show how this 
can be used to calculate A a and hence match with the resummed result. The 0(a s ) result, 
which comes from the ratio of the O{o i s ) three-parton and 0(a 2 ) two-parton cross sections, 
is only non-zero for P\ and P 2 - Furthermore, we have the relation Pi + P 2 = 1 + 0(a 2 ), 
which means that to this order we only need to calculate P2. We define the 0(a s ) term in 
the expansion of P2 to be R2 and obtain 

n( p , Jgjygg SXL drf ^ d<P CT dy 1%% dz^j 5 f 1 (x 1 )h(x 2 ) \M 3 \ 2 
my™t,R;Pt,r),<t>) - (27r) 4 ! ( r )f(r)\M P 

(47) 

in the pt-scheme J7|. The rescaled transverse momentum, y, energy fraction, z and azimuth 
around the jet axis, (p, are the most convenient variables to describe the two-parton jet (along 
with p t , r) and 4>) and are integrated over the region where the two partons are resolved as 
distinct subjets j7|. It is straightforward to perform these integrals numerically and hence 
compute R2. We can also easily classify the jets into 'quark', i.e. quark+gluon or anti- 
quark+gluon, 'gluon', i.e. gluon+gluon or quark+antiquark of the same flavour, or 'other' 
contributions. In all plots we show the sum of all three contributions, together with the 
separate quark and gluon contributions. 

To this order, the generating function is therefore given by 

$ jet (it, Q)=u(l- R 2 ) + u 2 R 2 + 0{a 2 s ). (48) 

In order to extract the leading-order probability of emitting a soft initial-state gluon into 
the jet, we take R2 and subtract from it the 0(a s ) expansion of the final-state contribution: 

P(y cut , R; p t , V, 0) = F q(Pu V, <f>)P q (Vcut, R, Pt, V, 0) + F g(Pt, V, ^)Pg(ycut, R; Pt, V, 0) 

= R 2 (y cm, R; Pt, v, 0) - A2(y C ut/R 2 ;pt,v,<fi), (49) 



where 



A 2 (F; p t , 7], 0) = F q (p t , 77, 0) A 2q (Y) + F g {p t , V , <f>)A 2g {Y), (50) 

A 2q (Y) = f Q 2 dQ' 2 Y q {Q,Q') + 0{a 2 s ) (51) 
= ^{C F \n 2 Y + 3C F lnF), (52) 

47T 

&2 9 (Y) = [ Q dQ' 2 (T g (Q,Q') + T f (Q')) + 0(a 2 s ) (53) 
JQl 

= ^(C A \n 2 Y + (3 \nY). (54) 

47T 
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The only remaining logarithmically-enhanced terms must come from initial-state radiation. 

According to the definition of A a , Q4"3]), the logarithmically-enhanced part of P is also 
equal to 

P a = a s A a \n(l/Y) + 0(a s ) } (55) 

where the neglected terms are pure 0(a s ) with no logarithmic enhancement. We can therefore 
equate P9| ) and fl5"5|) to logarithmic accuracy, to give 

_ (R 2 - A 2 ) a 

K ~ a s minY- (56) 

Inserting this into ([46]), we arrive at our central result for the generating function of jets in 
hadron collisions at scale Q = p t , 



- ui?° thers 4Vi£f hers . (57) 



We recall that $> a (u,Q) are computed in the NLLA as given in ( pl| , p^ ,p^-|35| ) and that 
d§ a (u, Q)/d\nY are computed in the LLA, by using the same formulas, but keeping only 
the logarithmically-enhanced terms in T a , (|24||25| ). 

Note that by computing R 2 using the full 2 — > 3 matrix elements, we are going beyond 
the NLLA. This is because we are guaranteed to compute the lowest order in a s contribution 
exactly, i.e. we have succeeded in matching the NLL approach with the fixed order calculation. 
We have therefore attained the accuracy aimed for in (|8|). The terms involving i?2 thers i n (0) 
are 0(a s ) terms that lie outside of the NLLA but are needed to formally ensure the matching. 
They are responsible for describing jets that are neither quark jets nor gluon jets and, as 
mentioned earlier, they are negligibly small. 



We finally note that when using the modified results of Section [2.2| , the modifications affect 
not only the resummed parts of fl5T|) , but also the matching term A 2 . Since this subtracts 
off the double-counting between the fixed-order term R 2 and the resummed results, and R 2 
already contains the correct threshold behaviour, A 2 must be calculated in the same way as 
the resummed results. Specifically, the inclusion of the Qo-suppressed terms in the evolution 
affects $ a and A 2 through equations (^T],[53|) and (|55|-|4Ti|) and the rescaling of Y affects all 
results except R 2 . 

In Figures |3] and ^ we show R 2 — A 2 for quark and gluon jets, without and with the 
threshold modifications. Since the latter are different for the different subjet rates, there are 
three solid curves corresponding to the one/two subjet rates (lower curve), three subjet rate 
(middle curve) and four subjet rate (upper curve). It can be clearly seen that for small Y, 
R 2 — A 2 consists of a NLL piece (~ a s InF) and a fixed order piece (~ a s ), i.e. 

R 2 - A 2 ~ a s ( InF + constant ). (58) 
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Figure 3: Comparison of (i?2 — A 2 ) g com- 
puted in the NLLA (dashed) with that com- 
puted after including threshold modifica- 
tions (solid). 
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Figure 4: Comparison of (i?2 — A 2 ) 9 com- 
puted in the NLLA (dashed) with that com- 
puted after including threshold modifica- 
tions (solid). 



Physically the logarithmically-rising term corresponds to initial-state radiation, while the 
constant term cannot uniquely be ascribed to the initial or final state. Note that the coefficient 
of the logarithm is unchanged by the threshold matching, while the constant term, which is 
small in the NLLA, becomes significant throughout the y cu t region we consider. We stress 
that this non-logarithmic part is included to ensure proper matching with the lowest order 
matrix element. We have chosen to multiply it by the LLA parton multiplication factor, 
<i$ a (w, Q)/d\n Y. In so doing we have introduced a tower of sub-leading terms that lie beyond 
our level of approximation. 

Looking at the structure of (0), it is tempting to think of the first term as arising solely 
from final-state radiation and the second solely from initial-state radiation. Indeed this is 
the case for logarithmically-enhanced terms, but constant terms like the one in (KSj) are not 
uniquely part of either contribution. The fact that i? 2 — A 2 is small at some y cn t value should 
not therefore be assumed to mean that initial-state effects are small there. 

Note that the slopes of Figures ^ and f|, i.e. the coefficients of the initial-state logarithms 
for quark and gluon jets, are rather similar. In fact in 0, it was mentioned that in the planar 
approximation (i.e. in the limit of large number of colours, iV c ) they would be the same. Their 
numerical difference (about 20%) is consistent with this statement. Moreover it was found 
that the slope was largely independent of their production mechanism, provided the rapidity 
of the other jet in the event was integrated out. If the recoiling jet's rapidity is constrained, 
the slope of i? 2 — A 2 becomes dependent on it. We return to this point later. 

It is also worth pointing out that the initial-state component is roughly proportional to R 2 , 
so reducing R can make the relative effect of initial-state radiation considerably smaller 0. 
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3 Numerical results 



We recall that our canonical jets are denned to be produced in pp collisions with p t = 60 GeV, 
T] = 0, R = 1 and \/s = 1800 GeV, and we use the CTEQ4M parton distribution functions [[L5| 
as implemented in PDFLIB (NSET = 34) 0. 



Canonical Jets 



Canonical Jets 
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Figure 5: One subjet rate: comparison 
of leading order (dotted), LLA (dashed), 
matched NLLA (solid) and matched NLLA 
with threshold modifications (green). 



Figure 6: Two subjet rate: comparison 
of leading order (dotted), LLA (dashed), 
matched NLLA (solid) and matched NLLA 
with threshold modifications (green). 




Figure 7: Three subjet rate: comparison of 
LLA (dashed), matched NLLA (solid) and 
matched NLLA with threshold modifications 
(green) . 



Figure 8: Four subjet rate: comparison of 
LLA (dashed), matched NLLA (solid) and 
matched NLLA with threshold modifications 
(green) . 



Figures ||-[8] show the individual subjet rates in canonical jets for up to and including four 
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subjets. Comparison is made between the leading order prediction, the prediction to leading 
logarithmic accuracy (LLA) and the fully matched next-to-leading logarithmic accuracy pre- 
diction (NLLA) with and without the improved treatment of the threshold region discussed 
in Sect. |2.2| . Firstly, we see that the fixed order results have only a small region of reliability, 
rapidly becoming unphysical for y cnt < 1CT 2 . The LLA result is physically behaved, but the 
higher subjet rates grow much too quickly away from threshold. The matched NLLA results 
on the other hand should be reliable for all y cut . They approximate the fixed-order results 
for Pi and P2 for large y cu t and remain well-behaved for small y cu t- We see that the modified 
threshold treatment makes very little difference for Pi, slightly more for P2 and more still for 
the higher subjet rates. This is an indication of the relative importance of neglected next-to- 
next-to-leading logarithmic effects and hence of the accuracy of the whole calculation. This 
difference would be reduced significantly by matching to the relevant n-parton tree-level ma- 
trix elements, probably to a similar level to P 2 , but to reduce the dependence further would 
require working to at least NNLLA. For all subsequent figures we use the threshold-improved 
results. 



Vs = 630 GeV and 1800 GeV Vs = 630 GeV and 1800 GeV 




Figure 9: y cut dependence of subjet rates in 
quark jets at 630 GeV (green) and 1800 GeV 
(black) and from final-state logarithms only 
(dashed). Also shown is the sum of all the 
rates. 



Figure 10: y cut dependence of subjet rates in 
gluon jets at 630 GeV (green) and 1800 GeV 
(black) and from final-state logarithms only 
(dashed). Also shown is the sum of all the 
rates. 



In Figure |9| we show the y cut dependence of the individual subjet rates in quark jets at the 
two different collider energies of 630 GeV and 1800 GeV. Figure |10] is a similar plot but for 



gluon jets and Figure |TT] is for all jets. Note the very weak dependence of the rates in quark 
and gluon jets on the centre-of-mass energy. This result supports the recent D0 analysis 
where it is assumed that jet observables do not depend upon centre-of-mass energy ||. One 
might be tempted to assume, seeing this result and the later ones, that because the properties 
of quark and gluon jets depend so little on how they were produced they are dominated by 
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Vs = 630 GeV and 1800 GeV 



Canonical Jets 
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Figure 11: y cut dependence of subjet rates in 
all jets at 630 GeV (green) and 1800 GeV 
(black) and from final-state logarithms only 
(dashed). Also shown is the sum of all the 
rates. 



Figure 12: Rates for iV subjets at 
y cut = 10~ 3 in quark (open circles) and 
gluon (solid circles) jets. The points in the 
N = 5 bin are for N > 5 subjets. 



final state effects. Comparison with the final-state-only curves in Figures |9] and [10] shows that 
this is not the case. A significant fraction of the resolved subjets come from non-final-state 
radiation but are nevertheless still universal to a very good approximation. The all-jets results 
do vary with centre-of-mass energy, because the quark-to-gluon mix is varying 



As discussed at the end of Section 2.4, the effect of initial state radiation cannot be 



inferred simply from the difference between the dashed and solid lines in Figures |9| and [10 



The difference represents the full effect of the R2 — A2 term in ([57]), after applying the 
threshold corrections. This contribution is made up of a logarithmic piece and an order a s 
piece (see (f)8|)). The logarithmic piece can be termed initial state radiation whilst the latter 
cannot. 

Figure [T2| shows the individual subjet rates at fixed y cut = 10" 3 for canonical quark and 
gluon jets. The point at iV = 5 is the inferred rate for 5 or more subjets. It has been computed 
using Pn>5 = 1 — Z)f=i wr th t ne first four subjet rates all shifted to have thresholds at Y 5 , 
so that Pn>5 is sensibly-behaved there. 



We turn now to the rapidity-dependence. Figures [13| and |14| show the y cut dependence of 
central (77 = 0) and forward (rj = 2) jets at p t = 60 GeV and y/s = 1800 GeV. Again they 
are almost indistinguishable, while the rates for all subjets, Figure [15|, do vary owing to the 
differing mix of quark and gluon jets. 

Figure [U] shows the y cut dependence of p t = 60 GeV and p t = 100 GeV jets at r\ = 
and yjs = 1800 GeV. This time the individual subjet rates for quark (Figure [17]) and gluon 
(Figure |TS|) jets do depend upon the jet p t , because the initial a s value is different. 
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Central and Forward Jets 



Central and Forward Jets 




Ycut 




Ycut 



Figure 13: Comparison of subjet rates for Figure 14: Comparison of subjet rates for 
rj = (black) and r\ = 2 (green) quark jets. rj = (black) and rj = 2 (green) gluon jets. 



Central and Forward Jets 60 GeV and 100 GeV Jets 




Figure 15: Comparison of subjet rates for 
r\ — (black) and r] = 2 (green) jets. 



Figure 16: Comparison of subjet rates 
for jets of p t = 60 GeV (black) and 
Pt = 100 GeV (green). 



Finally, we show in Figures |19] and ^0] that not only are the quark and gluon jet proper- 
ties independent of the collision energy, they are even independent of the collision type, by 
comparing canonical pp jets with those from 7r°7r° scattering at the same energy. To eliminate 
spurious differences due to different a s values in the different parton distribution functions, 
we choose, for these plots only, the GRV sets [JT7, I8fl , which use the same A value for both 
particle types. 

This independence of the properties of a given flavoured jet from the way it was produced is 
a highly non-trivial result. One would normally expect that the colour coherence of radiation 
from different emitters in an event would make the soft radiation dependent on the full details 
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60 GeV and 100 GcV Jets 60 GeV and 100 GeV Jets 




Figure 17: Comparison of subjet rates for 
quark jets of p t = 60 GeV (black) and 
p t = 100 GeV (green). 



Proton and Pion Quark Jets 




Ycut 



Figure 18: Comparison of subjet rates for 
gluon jets of p t = 60 GeV (black) and 
p t = 100 GeV (green). 



Proton and Pion Gluon Jets 




Ycut 



Figure 19: Comparison of subjet rates for 
quark jets in pp collisions (black) and in 7r°7r° 
collisions (green). 



Figure 20: Comparison of subjet rates for 
gluon jets in pp collisions (black) and in 7r°7r° 
collisions (green). 



of the hard scattering. As mentioned at the end of Section |27|, the amount of soft initial- 
state radiation into the jet is in fact largely independent of the scattering kinematics in the 
fully-inclusive case in which the recoiling jet is integrated out. However, putting requirements 
on additional jets breaks this inclusivity and the colour coherence changes the properties of 
the registered jet in response to the kinematics of the other jet. This is demonstrated in 
Figures ^T] and |22], where we show the quark and gluon subjet rates for jets at fixed rapidity 
in the forward region, 77 = +2, with the recoiling jet either unconstrained, or required to be 
at r] = +2 'same side' or r\ = — 2 'opposite side'. The 'same side/opposite side' ratio was 
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Figure 21: Comparison of sub jet rates for 
quark jets at 77 = 2 with the other jet in the 
event allowed to be at any rapidity (black), 
required to be at \t)\ = 2 on the same side 
(solid green) or opposite side (dashed green) 
of the event. 



Figure 22: Comparison of sub jet rates for 
gluon jets at rj = 2 with the other jet in the 
event allowed to be at any rapidity (black), 
required to be at \t)\ = 2 on the same side 
(solid green) or opposite side (dashed green) 
of the event. 



once thought to be a good way to separate quark and gluon jet properties since, for fixed jet 



kinematics, the quark-to-gluon mix is very different in the two samples. Figures [21] and |22 



show that this is not the the selection itself considerably biases the jet properties. 

Clearly the method of |J, which uses the centre-of-mass energy- dependence with fixed jet 
is superior. 

Although we have not performed a calculation for photoproduction [jp collisions), which 
has a slightly different structure owing to the direct photon contribution, we expect that 
the properties of inclusively-defined jets there would be similar to those in hadron collisions, 
which we have calculated. However, because of the colour coherence effect just mentioned, 
this statement is unlikely to be true once cuts are made on the other jets in the event. In 
particular, it is common to try to separate experimentally the direct and resolved photon 
events by imposing cuts on x 7 , the fraction of the photon's momentum that is reconstructed 
in the hardest two jets in the event. It seems likely that this cut will bias the jet properties 
sufficiently that our calculation cannot be used for a quantitative analysis. 



4 Conclusion 



We have calculated the sub jet rates in hadron collisions to next-to-leading accuracy in log- 
arithms of y cu t to all orders in a s , matched with the exact 0(a s ) result at large y cut . To 
this accuracy the contribution from initial-state radiation is essential. Nevertheless, it is still 
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possible to separate quark and gluon jets, and we find that their properties are almost com- 
pletely independent of their production mechanism, depending only on their p t , provided that 
they are defined fully inclusively As soon as additional cuts are placed on the event, the jet 
properties become dependent on the details of the hard scattering. 

Judging by comparisons of our results with and without threshold matching, which for- 
mally differ only by uncalculated next-to-next-to-leading logarithmic terms, we conclude that 
these terms are rather large, especially for the higher subjet rates. While this situation could 
certainly be improved by matching the n-subjet rate to the tree- level n+l-parton matrix ele- 
ment, further improvement will be extremely difficult. It is likely that, as in e + e~ annihilation, 
Pi will remain the best-calculated subjet rate. 
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